Effective mass and tricritical point for lattice fermions localized by a random mass 
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This is a numerical study of quasiparticle localization in symmetry class BD (realized, for example, 
in chiral p-wave superconductors) , by means of a staggered-fermion lattice model for two-dimensional 
Dirac fermions with a random mass. For sufficiently weak disorder, the system size dependence of 
the average (thermal) conductivity a is well described by an effective mass M^a, dependent on the 
first two moments of the random mass M{r). The effective mass vanishes linearly when the average 
mass M — >■ 0, reproducing the known insulator-insulator phase boundary with a scale invariant 
dimensionless conductivity ac — 1/n and critical exponent v = 1. For strong disorder a transition 
to a metallic phase appears, with larger ac but the same v. The intersection of the metal-insulator 
and insulator-insulator phase boundaries is identified as a repulsive tricritical point. 

PACS numbers: 72.15.Rn, 73.20.Jc, 74.25.fc, 74.78. Na 
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I. INTRODUCTION 
A. Description of the problem 

Superconductors with neither time-reversal symmetry 
nor spin-rotation symmetry (for example, having chiral 
p-wave pairing) still retain one fundamental symmetry: 
the charge-conjugation (or particle-hole) symmetry of 
the quasiparticle excitations. Because of this symmetry, 
quasiparticle localization in a disordered chiral p-wave 
superconductor is in a different universality class than in 
a normal metalji' The difference is particularly interesting 
in two dimensions, when the quantum Hall effect governs 
the transport properties. The electrical quantum Hall 
effect in a normal metal has the thermal quantum Hall 
effect as a superconducting analogue, with different 
scaling properties because of the particle-hole symmetry. 

The quasiparticles in a superconductor have electron 
and hole components ipejiph that are eigenstates, at ex- 
citation energy e, of the Bogoliubov-Dc Gcnnes equation 



Hq — Ep 
At 



A 

'Hq + Ep 



= e 



(1.1) 



In a chiral p-wave superconductor the order parameter 
A — ^{x{t'),Px — iPy} depends linearly on the momen- 
tum p = —ihd/dr, so the quadratic terms in the single- 
particle Hamiltonian Hq = /2m + U{r) may be ne- 
glected near p = 0. For a uniform order parameter 
xii") = AOj the quasiparticles are eigenstates of the Dirac 
Hamiltonian 



viPxUx +PyO-y) + v^M{r)az 



(1.2) 



with velocity v = xo and mass M — {U — Ep)/xo (dis- 
tinct from the electron mass m). The particle- hole sym- 
metry is expressed by 



Dirac ■ 



(1.3) 



Randomness in the electrostatic potential U (r) translates 
into randomness in the mass M{r) = M + SM{r) of 



the Dirac fermions. The sign of the average mass M 
determines the thermal Hall conductance^'"— which is 
zero for M > (strong pairing regime) and quantized 
at Go = TT^/c^T/G/i for M < (weak pairing regime). 

The thermal quantum Hall transition at Af = is 
analogous to the electrical quantum Hall transition at 
the center of a Landau level, but the scaling of the ther- 
mal conductivity a near the phase boundary is different 
from the scaling of the electrical conductivity because 
of the particle-hole symmetry. A further difference be- 
tween these two problems appear if the superconducting 
order parameter contains vorticesi^i^ A vortex of unit 
vorticity contains a Majorana bound state at e = 0, in 
the weak-pairing regime. ^'^ A sufficiently large density of 
Majorana bound states allows for extended states at the 
Fermi level (s ~ 0), with a thermal conductivity increas- 
ing cx InL with increasing system size LjS, This socalled 
thermal metal has no counterpart in the electronic quan- 
tum Hall effect. 

As indicated in Fig. [U there are two types of phase 
transitionSf^^ a metal-insulator (M-I) transition upon 
decreasing SM at constant M and an insulator-insulator 
(1-1) transition upon decreasing M through zero at con- 
stant (not too large) 6M. The I-I transition separates 
phases with a different value of the thermal Hall con- 
ductance, while the M-I transition separates the thermal 
metal from the thermal insulator. Only the I-I transi- 
tion remains if there are no vortices, or more generally, if 
there are no Malorana bound statesi^i^!^ In the nomen- 
clature of Ref. the symmetry class is called BD with 
Majorana bound states and D without. 



B. Objectives of this paper 

The primary purpose of our paper is to investi- 
gate, by numerical simulation, to what extent the 
scale dependence of localization by a random mass can 
be described in terms of an effective non-fluctuating 
mass: a{L,M,SM) ~ a{L, Mcs,0), for some function 
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FIG. 1: Phase diagram in symmetry class BD, calculated nu- 
merically from the lattice model of staggered fermions de- 
scribed in Sec. [Ill (A qualitatively similar phase diagram was 
calculated for a different model^^ in Refs. i and M ) The 
thermal conductivity decays exponentially oc e~^^^ in the 
localized phase and increases (x InL in the metallic phase. 
The thermal conductivity is scale invariant on the metal- 
insulator (M-I) phase boundary (red solid line), as well as 
on the insulator-insulator (I-I) phase boundary (blue dashed 
line). The M-I and I-I phase boundaries meet at the tricriti- 
cal point SM* . 



in (even) powers of M with the expansion of the scaling 
function (|1.4p in powers of L. This value for ly is aspect- 
ratio independent and agrees with the known result for 
the I-I transition.^ 

The description in terms of an effective mass breaks 
down for strong disorder. We find that the scaling func- 
tion at the M-I transition differs appreciably from ctq, 
with an aspect-ratio independent critical conductivity 
(Tc « 0.4 Gq. The critical exponent remains close to or 
equal to = 1 (in disagreement with earlier numerical 
simulations^*'). 

The secondary purpose of our paper is to establish the 
nature of the tricritical point SM* at which the two insu- 
lating phases and the metallic phase meet. From the scale 
dependence of a near this tricritical point, we conclude 
that it is a repulsive fixed point (in the sense that a scales 
with increasing L to larger values for 6M > SM* and to 
smaller values for 5M < 5M*). An attractive tricriti- 
cal point had been suggested as a possible scenario,— 
in combination with a repulsive critical point at some 
SM** < SM*. Our numerics does not support this sce- 
nario. 



Meg{M,SM). Because there is no other length scale in 
the problem at zero energy, cr(L, Meg, 0) can only depend 
on L and Mcs through the dimensionless combination 
LMcsv/h = L/S_. The effective-mass hypothesis thus 
implies one-parameter scahng: a{L,M,SM) = ao{L/^). 
Two further implications concern the critical conductiv- 
ity cTc (which is the scale invariant value of a on the phase 
boundary M = 0) and the critical exponent v (governing 
the divergence of the localization length ^ cx AI~'^). 

Both (Tc and ly follow directly from the effective mass 
hypothesis. By construction, the scaling function cto 
is the conductivity of ballistic massless Dirac fermions, 
which has been calculated in the context of graphene. For 
a system with dimensions L xW , and periodic boundary 
conditions in the transverse direction, it is given }yyl^^ 

J oo 

ao(L/0 = Go— cosh-2 ^[2i:nL/WY + {L/S,Y 

n— — oo 

^^Go- dq cosh-V?" + (i/C)"- (1-4) 

7^ Jo 

A scale invariant conductivity 

L ^ 

lim cro(i/C) = (Tc = Go— V cosh-^{2TTnL/W) 

^— )-oo [V ^ — ^ 

n— — OO 

(1.5) 

is reached for vanishing effective mass. In the limit of a 
large aspect ratio W j L ^ 1 we recover the known value 
<^c = Go /tt of the critical conductivity for a random mass 
with zero average."^'* The critical exponent v — \ follows 
by comparing the expansion of the conductivity 

(t(L, Af, SM) = (Tc + [Li/''M/(JAf)]2 -f- 0{Mf (1.6) 



C. Outline 

The outline of this paper is as follows. In the next Sec- 
tion we summarize the lattice fermion model that we use 
to simulate quasiparticle localization in symmetry class 
BD. We only give a brief description, referring to Refs.flTl 
and 18' for a more detailed presentation of the model. The 
scaling of the thermal conductivity and the localization 
length near the insulator-insulator and metal-insulator 
transitions are considered separately in Sees. IIIll and IIV) 
respectively. The tricritical point, at which the two phase 
boundaries meet, is studied in Sec. [V] We conclude in 
Sec. EH 



II. STAGGERED FERMION MODEL 

Earlier numerical investigationsSiiSii^ii^ of the class 
BD phase diagram were based on the Cho-Fisher net- 
work model. Here we use a staggered fermion model 
in the same symmetry class, originally developed in the 
context of lattice gauge theoryi^i^ and recently adapted 
to the study of transport properties in graphene?^ An 
attractive feature of the lattice model is that, by con- 
struction, it reduces to the Dirac Hamiltonian on length 
scales large compared to the lattice constant a. 

The model is defined on a square lattice in a strip 
geometry, extending in the longitudinal direction from 
X — QXox = L = Nxa and in the transverse direction 
from y = Qioy — W = NyU. We use periodic boundary 
conditions in the transverse direction. The transfer ma- 
trix T from X = to a; = L is derived in Ref . [TtI . and we 
refer to that paper for explicit formulas. The dispersion 
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relation of the staggered fermions, 



tan 



>(fcW2)+tan^(V/2)+(^)'=(£) 



ea 

(2.1) 

has a Dirac cone at wave vectors \k\a 1 which is 
gapped by a nonzero mass. Staggered fermions differ 
from Dirac fermions by the pole at the edge of Brillouin 
zone {\kx\ — 7r/a or \ky \ 7r/a), which is insensitive to 
the presence of a mass. We do not expect these large- 
wave number modes to affect the large-length scaling of 
the conductivity, because they preserve the electron-hole 
symmetry. 

The energy is fixed at e = (corresponding to the 
Fermi level for the superconducting quasiparticles) . The 
transfer matrix T is calculated recursively using a sta- 
ble QR decomposition algorithm.^ An alternative sta- 
bilization method^^ is used to recursively calculate the 
transmission matrix t. Both algorithms give consistent 
results, but the calculation of T is more accurate than 
that of t because it preserves the electron-hole symmetry 
irrespective of round-off errors. 

The random mass is introduced by randomly choos- 
ing values of M on each site uniformly in the interval 
(M - SM,M + SM). Variations of M(r) on the scale 
of the lattice constant introduce Majorana bound states, 
which place the model in the BD symmetry class In 
principle, it is possible to study also the class D phase 
diagram (without Majorana bound states), by choosing a 
random mass landscape that is smooth on the scale of a. 
Such a study was recently performed^^ using a different 
model,~ to demonstrate the absence of the M-I transi- 
tion in class D.^- Since here we wish to study both the 
I-I and M-I transitions, we do not take a smooth mass 
landscape. 



III. SCALING NEAR THE 
INSULATOR-INSULATOR TRANSITION 
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FIG. 2: Average conductivity a (with error bars indicating 
the statistical uncertainty) at fixed disorder strength SM — 
2.5h/va, as a function of system size L. The aspect ratio 
of the disordered strip is fixed at W/L — 5. Data sets at 
different values of M (listed in the figure in units of Ti/va) 
collapse upon rescaling by ^ onto a single curve (solid line), 
given by Eq. (|1.4|l in terms of an efi^ective mass M^s ~ h/v^. 
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FIG. 3: Plot of the average conductivity a versus {MLf, for 
fixed SM — 2.5 h/va and W/L = 3. The dashed line is a least- 
square fit through the data, consistent with critical exponent 
v = l. 



A. Scaling of the conductivity 

In Fig. [5] we show the average (thermal) conductiv- 
ity a = {L/W){Titt^) (averaged over some 10^ disorder 
realizations) as a function of L for a fixed 6M in the 
localized phase. Data sets with different M collapse on 
a single curve upon rescaling with ^. (In the logarith- 
mic plot this rescaling amounts simply to a horizontal 
displacement of the entire data set.) The scaling curve 
(solid line in Fig. is the effective mass conductivity 
(|1.4p . with Mctf — h/v^. Fig. [3] shows the linear scaling 
of (J with [MLY for small A/, as expected from Eq. (II. 6p 
with V — \. 

We have studied the aspect ratio dependence of the 
critical conductivity CTc. As illustrated in Fig.|4l the con- 
vergence for W/ L — oo is to the value CTc = I/tt expected 
from Eq. (II. 4p . The conductivity of ballistic massless 
Dirac fermions also has an aspect ratio dependenceflS 



given by Eq. (jl.5p (for periodic boundary conditions). 
The comparison in Fig. 0] of CTc with Eq. (|1.5p shows that 
(Tc has a somewhat weaker aspect ratio dependence. 

B. Scaling of the Lyapunov exponent 

The transfer matrix T provides an independent probe 
of the critical scaling through the Lyapunov exponents. 
The transfer matrix product 7~T^ has eigenvalues e^^" 
with < /^i < /i2 < • • • • The n-th Lyapunov exponent 
a„ is defined by 

an = lim (3.1) 

The dimensionless product Wai = A is the inverse of 
the MacKinnon-Kramer parameter.?^ We obtain ai by 
increasing L at constant W until convergence is reached 
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FIG. 4: Dependence on the aspect ratio W/L of the critical 
conductivity at the insulator-insulator transition (M = 0). 
The disorder strength is fixed at 5M = 2.^'h/va. The dashed 
curve is the aspect ratio dependence of the conductivity of 
ballistic massless Dirac fermions [Eq. p.Sp ]. 
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FIG. 5: Plot of A = Wa\ (with ai the first Lyapunov expo- 
nent) as a function of M near the insulator-insulator transi- 
tion, for fixed 5M = 2.5hv/a and different values of W . The 
dashed lines are a fit to Eq. (13. 2p . 
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FIG. 6: Average conductivity a at fixed average mass M = 
0.032 ?i/«a, as a function of system size L. (The two panels 
show the same data on a different scale.) The aspect ratio of 
the disordered strip is fixed at W/L = 5. Data sets at different 
values of SM (listed in the figure in units of h/va) collapse 
upon rescaling by ^ onto a pair of curves in the metallic and 
insulating regimes. The metal-insulator transition has a scale 
invariant conductivity a'^, larger than the value Go/vr which 
follows from the effective mass scaling (dashed curve in the 
lower panel). The upper panel shows that the conductivity in 
the metallic regime follows the logarithmic scaling (|4.1|l . 



(typically for L/W ~ 10^). The large-L limit is self- 
averaging, but some improvement in statistical accuracy 
is reached by averaging over a small number (10-20) of 
disorder realizations. 

We seek the coefficients in the scaling expansion 

A = Ac + ciW^^^iM - Mc) + 0{M - M^f, (3.2) 

for fixed 5M. The fit in Fig. Ogives A^ = 0.03, u = 1.05, 
Mc ~ 7 • 10~^, consistent with the expected values^ Ac = 
0, = 1, Mc = 0. 



IV. SCALING NEAR THE METAL-INSULATOR 
TRANSITION 

A. Scaling of the conductivity 

To investigate the scaling near the metal-insulator 
transition, we increase 5M at constant M . Results for 
the conductivity are shown in Fig. [51 In the metallic 



regime 5M > SM^ the conductivity increases logarithmi- 
cally with system size L, in accord with the theoretical 
prediction)ii^ 

(t/Go = — hiL + constant. (4.1) 

(See the dashed line in Fig.[6l upper panel.) 

In the insulating regime SM < 5Mc the conductivity 
decays exponentially with system size, while it is scale 
independent at the critical point 6M = SMc- Data sets 
for different 5M collapse onto a single function of i/^, 
but this function is different from the effective mass scal- 
ing (Jq{L/S) of Eq. (|1.4p . (See the dashed curve in Fig.lH 
lower panel.) This indicates that the effective mass de- 
scription, which applies well near the insulator-insulator 
transition, breaks down at large disorder strengths near 
the metal-insulator transition. The two transitions there- 
fore have a different scaling behavior, and can have differ- 
ent values of critical conductivity and critical exponent 
(which we denote by a'^ and v'). 

Indeed, the critical conductivity a'^ = 0.41 Gq is sig- 
nificantly larger than the ballistic value Gq/tt = 0.32 Gq. 
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FIG. 7: Plot of the average conductivity cr as a function of SM 
near the metal-insulator transition, for fixed M — 0.032 h/va. 
The length L is varied at fixed aspect ratio W/L = 3. The 
dashed curves are a fit to Eq. (14. 2|) . 
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FIG. 8: Plot of A = Wai (with ai the first Lyapunov expo- 
nent) as a function of SM near the metal-insulator transition, 
for fixed M = 0.032 hv/ a and different values of W. The 
dashed curves are a fit to Eq. (|4.3p . 

Unlike at the insulator-insulator transition, we found no 
significant aspect-ratio dependence in the value of cr^ . To 
obtain the critical exponent i^' we follow Ref. [25| and fit 
the conductivity near the critical point including terms 
of second order in SM — SAI^. 

(T = cr^ + ciL^/"' [SM - SMa + C2{5M - SM^)^] 

+ c^L^'''\5M -SKhf. (4.2) 

Results are shown in Fig. H with v' = 1.02 ± 0.06. The 
quality of the multi-parameter fit is assured by a reduced 
chi-squared value close to unity (x^ = 0.95). Within 
error bars, this value of the critical exponent is the same 
as the value v = 1 for the insulator-insulator transition. 



B. Scaling of the Lyapunov exponent 

As an independent measurement of v' , we have inves- 
tigated the finite-size scaling of the first Lyapunov expo- 
nent. Results are shown in Fig. [51 Within the framework 



FIG. 9: Conductivity cr as a function of SM on the critical 
line M = 0, for different values of L at fixed aspect ratio 
W/L = 3. (The dotted lines through data points are guides 
to the eye.) The tricritical point SM* is indicated, as well as 
the scale invariant large-L limit Gc for SM < SM* . 

of single-parameter scaling, the value of v' should be the 
same for a and A, but the other coefficients in the scaling 
law may differ, 

A = Ac + c'^L^^"' [SM - SM'^ + c'^iSM - SM'J^] 

+ c'^L^''''{SM ~SM'^^. (4.3) 

Results are shown in Fig. [HI with v' = 1.06 ± 0.05. The 
chi-squared value for this fit is relatively large, = 5.0, 
but the value of v' is consistent with that obtained from 
the conductivity (Fig. [7]). 

V. TRICRITICAL POINT 

As indicated in the phase diagram of Fig. [TJ the tri- 
critical point at AI = 0, SM = SM* is the point at which 
the insulating phases at the two sides of the I-I transi- 
tion meet the metallic phase. We have searched for this 
tricritical point by calculating the scale dependence of 
the conductivity a on the line M = Q for different SAL 
Results are shown in Fig. [5] 

The calculated scale dependence is consistent with the 
identification of the point SM* = 3.44 h/va as a repulsive 
fixed point. The conductivity increases with increasing L 
for SAI > SAI* , while for SM < SM* it decreases towards 
the scale invariant large-L limit CTg. 

VI. DISCUSSION 

We have studied quasiparticle localization in symme- 
try class BD, by means of a lattice fermion model. The 
thermal quantum Hall effecli^^— in a chiral p-wave super- 
conductor at weak disorder is in this universality class, 
as is the phase transition to a thermal metal^ at strong 
disorder. For weak disorder our lattice model can also 
be used to describe the localization of Dirac fermions in 
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graphene with a random gap^^i^i^l (with a the electri- 
cal, rather than thermal, conductivity and Gq = 4,e^/h 
the electrical conductance quantum) . The metallic phase 
at strong disorder requires Majorana bound statesj^i^ 
which do not exist in graphene (symmetry class D rather 
than BD). Wc therefore expect the scaling analysis in 
Sec. mil at the insulator-insulator (I-I) transition to be 
applicable to chiral p-wave superconductors as well as 
to graphene, while the scaling analysis of Sec. |lV]at the 
metal-insulator (M-I) transition applies only in the con- 
text of superconductivity. 

Our lattice fermion model is different from the network 
modeUi used in previous investigations ^ i is , i6 
falls in the same universality class so we expect the same 
critical conductivity and critical exponent. For the I- 
I transition analytical calculations^'^* give CTc = Gq/tt 
and J/ = 1, in agreement with our numerics. There are 
no analytical results for the M-I transition. We find a 
shghtly larger critical conductivity (cr^ — 0.4 Go), which 
has the qualitatively more significant consequence that 
the effective mass scaling which we have demonstrated 
at the I-I transition breaks down at the M-I transition 
(compare Figs. [H and [51 lower panel). 

We conclude from our numerics that the critical expo- 
nents V at the I-I transition and v' at the M-I transition 
are both equal to unity within a 5% error margin, which is 
significantly smaller than the result — v' — 1.4 ± 0.2 of 
an earlier numerical investigation^^ It is much more diffi- 
cult to calculate accurately the critical exponent than the 
conductivity, so we do not attach too much significance 
to this disagreement. (Athough we find it puzzling that 



the value for i/ obtained in Ref. [13 lies so much higher 
than the analytical prediction v = 1.) The logarithmic 
scaling (|4.ip of the conductivity in the thermal metal 
phase, predicted analytically?^ is nicely reproduced by 
our numerics (Fig.|6l upper panel). 

The nature of the tricritical point has been much de- 
bated in the literatureJ^ii^ Our numerics indicates that 
this is a repulsive critical point (Fig. [5]). This finding 
lends support to the simplest scaling fiow along the I-I 
phase boundary,i4 towards the free-fermion fixed point 
at M = 0,SM = 0. 

In conclusion, we hope that this investigation brings 
us further to a complete understanding of the phase dia- 
gram and scaling properties of the thermal quantum Hall 
effect. We now have two efficient numerical models in the 
BD universality class, the Cho-Fisher network modelii 
studied previously and the lattice fermion modeli^ stud- 
ied here. There is a consensus on the scaling at weak 
disorder, although some disagreement on the scaling at 
strong disorder remains to be resolved. 
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